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ABSTRACT 

We revisit the flip-flop instability of two-dimensional planar accretion using high-fidelity numerical 
simulations. By starting from an initially steady-state axisymmetric solution, we are able to follow the 
growth of this overstability from small amplitudes. In the small-amplitude limit, before any transient 
accretion disk is formed, the oscillation period of the accretion shock is comparable to the Keplerian 
period at the Hoyle-Lyttleton accretion radius {Ra), independent of the size of the accreting object. 
The growth rate of the overstability increases dramatically with decreasing size of the accretor, but is 
relatively insensitive to the upstream Mach number of the flow. We confirm that the flip-fiop does not 
require any gradient in the upstream flow. Indeed, a small density gradient as used in the discovery 
simulations has virtually no influence on the growth rate of the overstability. The ratio of specific 
heats does infiuence the overstability, with smaller 7 leading to faster growth of the instability. For 
a relatively large accretor (a radius of 0.037 ii!o) planar accretion is unstable for 7 = 4/3, but stable 
for 7 > 1.6. Planar accretion is unstable even for 7 — 5/3 provided the accretor has a radius of 
< 0.0025 i?a. We also confirm that when the accretor is sufficiently small, the secular evolution is 
described by sudden jumps between states with counter-rotating quasi-Keplerian accretion disks. 
Subject headings: accretion — hydrodynamics — shock waves — turbulence 



1. INTRODUCTION 

The instability of hydrodynamic accretion from a mov- 
ing medium and its relevance to the temporal behav- 
ior of accreting x-ray sources has been debated since 
the discovery of a fiip-fiop instability in two-dimensional 
hydr odynamic simulations o f an accreting compact ob- 
ject dFryxell fc Taam| 11988] hereafter FT). This non- 
axisymmetric instability has been subsequently seen in 
numerous hydrodynamic simulations of t wo-dimensional 
planar accretion ([M atsuda et al. 1991: Zarinelli et al. 
19951 iBen ensohn et al..l997;|Shimaiet al..l998. : Pogorelov 
et al. 2000|, but has not been un ambiguously i dentified 



m three-dimensional simulations (Ruffert 19991. Given 



that accretion-powered X-ray sources in Be-type binary 
systems are expected to accrete from a planar disk ex- 
pelled by the companion Be star, attention has focussed 
on the relationship of t he flip-flop instabil ity to these 
Be-type X-ray pulsars. 'Ta am et al.| ( 1988 1 interpreted 
the recurrent flares in EXO 2030-1-375 in terms of the 
episodic accretion resulting from the flip-flop instability. 
Benensohn et al. ( 1997j ) argued that the variable accre- 
tion torques associated with the flip-flop instability may 
explain the erratic spin behavior in these systems. 

The flrst quantitati ve description of hyd r odynam ic ac- 
cretion was given by 'Hoyle & Lyttleton' f 1939a'), who 
were attempting to explain climate changes on Earth in 
terms of an increase in solar output due to the addition of 
accretion energy as the Sun passed through an interstel- 
lar cloud. Shortly thereafter the same authors proposed 
hydrodynamic accretion from interstellar g solu- 
tion to the paradox of the extremely short lifetimes of 
massive stars compared to geologic timescales, which at 
the time was considered "an exceedingly embarrassing 
consequence of the presen t theory of stellar evolution" 
( |Hoyle fc Lyttleton|[l939bl ). 



The Hoyle-Lyttleton accretion model is parameterized 
by an object of mass M moving at a speed Voo through 
a uniform medium of densit y pop. Using ballistic orbits 
(neglecting pressure effects), Hoyle & Lyttleton (1939a I 
deflned an accretion radius given by 



Ra 



2GM 

1/2 ' 



(1) 



such that gas approaching the star with an impact pa- 
rameter less than Ra would collide on an accretion line 
behind the star and would be left with insufficient kinetic 
energy to escape the gravitational potential of the star. 
The mass accretion rate is expected to be comparable to 
the mass flux through a circle of radius Ra far upstream: 
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(2) 

A thorough review of Hoyle-Lyttleton accretion theory 
and s ubsequent numerical studies is given by [Edgar 
(j2004|. We provide only a brief synopsis here. 

This basic model of hydro dynamic acc retion was flrst 
investigated numerically by Hunt (19711, who used 2D 
axisymmetric steady-state solutions to quantify the mass 
accretion rate. For an upstream Mach number of 2.4 
he found an accretion rate of O.SSMhl- Subsequent 
time-dependent 2D axisymmetric simulations have, for 
the most part, sho wn that axisymmetric hydrodynamic 
accre ti on is stable (iShima et al. 1985 Pogorelov et al. 
20001. Fryxell et al. 1 1987p found steady axisymmetric 
accretion flows provided that the inner b o undar y con- 
dition was totally absorbing. Koide et al. ( 1991 1 found 
steady solutions for low Mach number flows, but they 
found the accretion flow to be unsteady if the upstream 
Mach number was greater than 2. These latter results, 
however, are likely plagued by poor resolution near the 
accreting object. The authors attributed the unsteady 
flow to gas flowing towards the accretor within the ac- 
cretion column, but slightly missing the accretor. One 
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might expect the gas in the accretion column to be com- 
pletely absorbed if a higher grid resolution near the ac- 
cretor fully resolved the flow. 

This picture of steady, axisymmetric accretion is al- 
tered in the context of an x-ray binary system in which 
a compact object is accreting from the stellar wind of a 
companion star. In this case gradients in both density 
and velocity of the flow past the compact star will break 
the axisymmetry, and in particular may lead to the accre- 



tion of an gular momentum as w ell as mass ( Illarionov fc 
SunyaevI 1975,: Shapiro & Lightman 19761. Matsuda et 



aT (19871 used 2D hydrodynamic sunuiations to study 



tEe dynamics of such wind accretion in a close binary 
system. They presented some models for which the flow 
near the accreting star sporadically changed its direction 
of rotation - a phenomenon the authors referred to as a 
'flip-flop flow.' 

Shortly thereafter, numerical studies of the more ideal- 
ized problem of Hoyle-Lyttleton accretion in 2D planar 
flow (p'T) found a similar behavior, with the accretion 
wake flipping from side to side. The result of this flip- 
flopping is the formation of temporary accretion disks 
spinning in alternating directions, with a burst of mass 
accretion when the direction of rotation is flipped. This 
phenomenon has been seen in many subsequent simula- 
tions, and has come to be known as the 'flip-flop' insta- 
bility. 

While Hoyle-Lyttleton accretion is a very idealized 
problem, it is important to note that the flip-flop insta- 
bility was flrst seen in more realistic si mulations of gas 



flow within a binary star system (Matsuda et al. 19871 



and in subsequent simulations that incorporated several 
key physical pr ocesses including x-r ay heating and ra- 



diative cooling (Blondin et al. 19901. Although in the 



latter case it was found that strong Compton heating 
due to a high x-ray luminosity of the accreting star sup- 
pressed the flip-flop instability. However, both of theses 
studies of accretion in binary systems were restricted to 
2D. The increase in density and pressure as the accre- 
tion flow approaches the compact star will be different 
in three dimensions compared to two, raising the possi- 
bility that the stability of 3D accretion may differ from 
these 2D results. The only published simulations of 3D 
Hoyle-Lyttleton accretion with a reasonably small accre- 
tor were a series of studies using a ne sted cartesian grid 
( [Ruffert fc Arnet t||1994| | Ruffert||l999[ ). The accretion in 
these 3D studies was found to be unstable, but the evolu- 
tion did not resemble the flip-flop behavior found in 2D 
planar flow. In particular, the amplitude of variations 
remained considerably smaller than that found in 2D, 
and no evidence of steady quasi-Keplcrian disk flow was 
ev er seen. Nonethele s s, som e of the 3D models reported 



Ruffert fc Arnett (19941 did exhibit a quasi-periodic 
swinging of the accretion shock leading to brief epochs of 
strong rotational flow near the accretor. This behavior 
was particularly evident for the models with the smallest 
accretors. Although these 3D runs used nested grids to 
improve spatial resolution near the accreting star, each 
nested grid was only 32^ zones. This resolution is rela- 
tively low compared to that used in many of the 2D simu- 
lations exhibiting the flip-flop, and may have contributed 
to the lack of detection of the flip-flop instability in 3D. 

Independent of the question of whether adiabatic, 2D 
planar flow is an appropriate model for astrophysical 



systems, the flip-flop instability itself has drawn signif- 
icant scrutiny and criticism over the past two decades. 
Many examples of 2D planar accretion simulations have 
been presente d in the literature, wit h mixed results as 



tabulated by Foglizzo et al. (20051. Several authors 
have pointed out various shortcommgs of existing 2D 
planar Hoyle-Lyttleton simulations. These include, but 
are not necessar ily limited to, poor numerical resolution 
([Zarinelli et al. 1995), inaccurate algorithms (Pogorelov 
e^^^Ofl |, inappropriate boundary conditions upstream 
and at the surface of the accretor ( jZarinelli et al.,,1995| 
Benensohn et al. 1997|), la ck of conservation of angular 
momentum ( |Shima et al. 19981, and an unrealistically 
large accretor. T'hese criticisms continue to cast doubt 
and confusion over the existence of the flip-flop instabil- 
ity. 

In this paper we present an in-depth study intended to 
provide a better understanding of the properties of the 
flip-flop instability. We take care to address the many 
issues raised about previous numerical simulations and 
we use the increased availability of computing power to 
more fully explore parameter space. We describe our 
numerical method in Section 2, building from the orig- 
inal flip-flop simulations of |FT[ We present our results 
in Section 3, where we show the flip-flop instability is a 
true overstability of the accretion wake, and provide a 
summary of our conclusions in Section 4. 

2. COMPUTATIONAL METHOD 

We use the time-dependent hydrodynamics code VH-1 
to solve the Euler equations for an ideal gas character- 
ized by an adiabatic index 7 in 2D polar coordinates 
with a point source of gravity at the coordinate origin. 

Both codes 



This code is similar to the code used by |FT 
are based on the Pi ece-wise Parabolic Method ( [Colella 
fc Woodward 19841, but the code used here uses the 
Lagrange-remap formulation while the code used by |FT| 
used the direct-Eulerian formulation. None-the-less, we 
expect very similar results from these two codes. 

To this end, we begin by presenting results for a simu- 
lation set up to mimic a simulation in |FT[ using the same 
grid, boundary conditions and simulation parameters. 
Specifically, we use a numerical grid with 100 radial zones 
and 200 azimuthal zones. The inner radius correspond- 
ing to the surface of the accreting star is Rg — 0.037 i?a, 
the first zone has a width of AR = 0.25 Rs, and the sub- 
sequent radial zones increase in width by a factor of 1.03. 
This puts the outer radius of the computational grid at 
6.1 i?a. The flow at the upstream boundary is given a 
uniform velocity with a Mach number of 4. The gradients 
for all flow variables at the outer boundary were set to 
zero for the downstream half of the grid, with the further 
condition that the radial velocity remained positive. We 
imposed an absorbing boundary condition at the surface 
of the accreting star by setting the density and pressur e 
to very small values inside of this boundary ( Hunt|197r i. 
We note that the flip-flop instability has been shown to 
be relatively insensitive to the speciflc boundary condi 
tions imposed on the accreting surface ( Taam fc Fryxell 
1989 



J jZarincUi et al. 1995). We imposed a transverse 



density gradient at the upstream boundary of the form 

(1 + ^ — e), (3) 



with e = 0.005. 
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Fig. 1. — The Holyle-Lyttelton accretion model with uniform 
flow at Mach 4 entering from the right and passing a point source 
of gravity, which is located at the apex of the shock cone marked 
by the heavy solid line. The gray lines are streamlines as computed 
in the steady-state initial conditions described in Section 2.3. The 
outermost streamlines correspond to an upstream impact param- 
eter of slightly greater than /?„ . The parameters of the model 
shown here are 7 = 4/3 and Rs = 0.0125 -Ra 




30 40 50 
time (RaA/„) 

Fig. 2. — The time evolution of the specific angular momen- 
tum of the accreted gas in our initial simulation with the standard 
VH-1 hydrodynamics code is qualitatively similar to the results of 
FT. The flip-flop instability is characterized by episodes of disk 
accretion with alternating directions of rotation. 



The results of this simulation, shown in Figure [2] as a 
trace of the specific angular momentum of the accreted 
gas, is qualitatively similar to the corresponding simula- 
tion presented in FT. (Note that we normalize the spe- 
cific angular momentum to that of a circular orbit at the 
radius of the accretor, RgVc — y/GMRg, so that the satu- 
ration value when a temporary accretion disk has formed 
is independent of the accretor radius.) In particular, the 
accretion is marked by phases of disk accretion with rel- 
atively constant specific angular momentum correspond- 
ing roughly to that of a circular orbit at the surface of 
the accretor. These disk phases occur for both directions 
of rotation, with sometimes very abrupt transitions be- 
tween the two. 

Given the dramatic increase in available computing 
power in the two decades since the discovery of the fiip- 



fiop instability, we are at liberty to run a large number 
of simulations to refine our computational model, to run 
at higher spatial resolution, and to explore the parame- 
ter space more completely. In the following subsections 
we delineate the changes made to this initial numerical 
model in order to arrive at what we consider to be high- 
fidelity hydrodynamic simulations demonstrating the ex- 
ponential growth of the fiip-fiop instability. 

2.1. Hydrodynamic Algorithm 

A serious shortcoming of our initial model was the pres- 
ence of density inhomogeneities in the upstream fiow due 
to numerical noise. The source of these numerical fiuc- 
tuations was traced to the conservation of total energy 
at the expense of ensuring constant entropy in the up- 
stream supersonic flow traveling at various angles with 
respect to the computational grid. To reduce errors in in- 
ternal energy created by subtracting kinetic energy from 
the conserved total energy, we modified the algorithm 
to remap internal energy rather than total energy when 
the Mach number of the flow exceeds 3.0. This modifi- 
cation still provides for accurate shock jump conditions 
because the original conservation of total energy is used 
in zones containing a shock, but it now allows one to 
model uniform supersonic (Mach = 4) planar fiow across 
a cylindrical polar grid with minimal numerical artifacts. 

2.2. Upstream Boundary Conditions 

Assuming a uniform planar fiow at the upstream 
boundary neglects any gravitational effects of the fiow 
from infinity to that boundary. As a result, the flow near 
the accretor becomes dependent on the chose n radius of 



the outer boundary (Benensohn et al. 1997). A more 



realistic boundary condition is to assume a ballistic tra- 



jectory from infinity to the outer boundary ( Matsuda et 



al. 1991 Benensohn et al. 1997 Shima et al. 



imr. 



result ing flow velocity at a given po int (r, 0) m space is 



then ( Bisnovatyi-Kogan et al.|[T979 ) 



Vr 



c 



(4) 



(5) 



where ( is the impact parameter at infinity of the stream 
line passing through that point. 



r2 sin^ - 2ri?a(l + cos < 



(6) 



One can then find the local gas density using mass 
conservation, equating the mass fiux through a region 
bounded by ^ and C -I- c?C infinity to the mass fiux 
through a region bounded by r and r + dr for a point 
along the same stream line. For the 2D planar geometry 
used in this paper the result is 



Po 



c 



2C 



(7) 



Once the density is known, the pressure can be found 
by assuming the fiow is isentropic. The fiow variables 
upstream of the accretion shock in our simulations re- 
main within 0.5 percent of this analytic solution for an 
upstream Mach number of 4. 
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2.3. Steady Initial Conditions 

Although often not exphcitly stated, previous simula- 
tions were typically started from an initial state of per- 
fectly uniform planar flow. The early evolution is then 
dominated by the formation of a tail shock along the 
symmetry axis behind the accretor and the subsequent 
expansion of the region of shocked gas. This highly vari- 
able initial flow masks the origin of any possible instabil- 
ity, and indeed makes it impossible to assess whether the 
observed behavior is a true instability in the sense that 
an initial equilibrium state will tend to evolve away from 
equilibrium. 

We address this problem by creating a steady-state 
solution for initial conditions. We evolve our model on 
only half of the grid, enforcing reflection symmetry about 
the line of centers. These 'half simulations are initialized 
with the ballistic flow model and allowed to evolve for 
several flow crossing times, during which the accretion 
shock settles into a steady state. These relaxed solutions, 
an example of which is shown in Figure [l] are then used 
to initialize the full 2D simulations. 

2.4. Downstream Boundary Conditions 

The downstream boundary condition in previous nu- 
merical simulations has traditionally been a condition of 
zero gradients, allowing flow to smoothly move off the 
grid. This is a reasonable approach if the flow is moving 
supersonically off the grid, but can lead to spurious re- 
ffections if the flow across the boundary is subsonic. To 
investigate the downstream conditions in the accretion 
wake we examined the Mach number of the flow both 
in our steady-state initial conditions and in the time- 
dependent unstable flow. In the steady-state flow the 
radial velocity on the symmetry axis does not become 
positive until ^ [width = 3.5m] Ra, so we focus our at- 
tention to flow slightly off this axis. In Figure |3] we plot 
the Mach number of the radial flow as a function of dis- 
tance from the accretor for angles of 7 and 14 degrees 
from the symmetry axis. For reference, the shock is at an 
angle of about 21 degrees far downstream. At 7 degrees 
from the axis, the flow exceeds Mach 1 only after trav- 
eling ~ 8i?a, and even then it barely rises above unity. 
We also show the minimum Mach number inside of the 
accretion wake for a simulation subject to the flip-flop in- 
stability. The unsteady flow falls in an envelope defined 
by the 7 and 14 degree traces. S etting the outer radius t o 



greater than 10 as done by Benensohn et al 
ensures that the gas will propagate supersonical 



(|1997f, 
y oti of 



the grid under most circumstances, thereby minimizing 
wave reflection from the boundary. 

2.5. Conserve Angular Momentum 



Shima et al. (19981 demonstrate that the dynamics of 



the temporary accretion disks formed during the flip-flop 
instability are less dependent on the underlying computa- 
tional mesh if the hydrodynamic algorithm strictly con- 
serves angular momentum rather than linear momentum. 
This result is not unexpected given that during the disk 
epochs mass accretion is driven by spiral shocks within a 
roughly Keplerian flow, and any mass accretion is closely 
tied to the transport of angular momentum. To this end 
we made a minor change to the PPMLR algorithm used 
by VH-1 in order to explicitly conserve angular momen- 
tum. This code modiflcation has been used to quantify 




Radius 

Fig. 3. — Mach number of the outward radial flow inside the 
accretion wake as a function of distance from the accretor (in units 
of Ra)- The two dashed Hues are from a steady-state solution, 
tracing the Mach number along angles of 7 and 14 degrees from 
the axis. The solid lines show the minimum Mach number inside 
the wake at different times in a simulation of unstable accretion 
flow. 



disk accretion driven by spiral shocks (Blondin 20001 



measuring effective values of the disk a parameter as 
small as 10~^. 

2.6. Spatial Resolution 

The computational grid used by |FT| resulted in zones 
with extreme aspect ratios (8:1) near the surface of the 
accretor, while other r esearchers used a sin gle zone to 
represent the accretor (Matsuda et al. 19911. More re- 



cent work has attempted to improve this shortcoming by 
running identical si mulations with increasing numerical 
resolution (B enensohn et al.||1997 l or by applying adap- 
tive mesh refinem ent to more accurate ly resolve the flow 
near the accretor ( Zarinelli et al.|1995 ). We chose to use 
a grid similar to that of |i:''T[ namely a uniform grid in 
(j) and a logarithmically increasing g rid in r. However, 



we followed the approach of Benenso hn et al.| (1997l and 
chose parameters that kept the aspect ratio of all zones 
roughly square. This means both that Ar/r is constant 
throughout the grid, and also that Ar is much smaller 
than the radius of the accretor given that we use a mini- 
mum of 200 zones in the angular direction and we require 
that Ar be comparable to rAi/i. 

3. RESULTS 

For our standard model (A in Table 1) we chose pa- 
rameters comparable to that of FT , namely an upstream 
Mach number of four, R^ = 0.037 i?^, and 7 = 4/3. We 
did not impose any asymmetry in the upstream flow. 
This simulation does exhibit the flip-flop instability (Fig- 
ure H , however the oscillation of the accretion wake is 
less abrupt than found in the original simulations of |FT[ 
Much of the evolution is characterized by the accretion 
wake swinging back and forth, but remaining attached to 
the accretor. Only for brief epochs does the shock wrap 
around the leading edge of the accretor and 'detach', re- 
sulting in a temporary accretion disk as seen in the last 
frame of Figure [4] Moreover, the flow never exhibits an 
abrupt transition from a disk of one rotation directly to a 
counter-rotating disk. Instead, the disk disappears and 
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the accretion wake swings back and forth one or more 
times before forming another disk. 

The mass accretion rate and the specific angular mo- 
mentum of the accreted gas are shown in Figure [5] The 
mass accretion rate is normalized to the 2D planar equiv- 
alent of equation (|2| : 



M2D ^ 2 RaPooVoc 



(8) 



During the initial evolution, before the oscillations of the 
accretion shock grow to the point of driving the shock 
around the leading edge of the accretor, the mass accre- 
tion rate remains at a roughly constant value of 0.9 M2D- 
This is i n good agreement with the original calculation 
of , Hunt] ( |1971 1 for axisymmetric Hoyle-Lyttleton accre- 
tion as well as subseque nt simulations of planar accretion 
flow ( [Koide et al.||1991 ) . During the brief epochs of disk 
accretion, the mass accretion rate is both substantially 
lower and highly variable. The termination of each disk 
phase is marked by a sudden spike in the mass accretion 
rate as the bulk of the mass that had been orbiting in a 
quasi-Keplerian accretion disk is rapidly accreted on the 
relatively short free-fall timescale. 

To quantify the growth of the instability we fit the 
early evolution (before the oscillations reach a maximum 
amplitude) with an exponentially growing sinusoid: 



j(t) = jo d^^*- cos{uj^t) 



(9) 



The fitted function shown in Figure [5] corresponds to 
an oscillation frequency of uoi = 0.68 and an exponential 
growth rate of uj^. = 0.07. This corresponds to the ampli- 
tude approximately doubling each period of oscillation. 
Note that the oscillation period {2-K/ijJi = 9.2) is com- 
parable to the Keplerian period at the accretion radius 
{2tt^/2 = 8.9). Eventually some of the peaks in the trace 
of j{t) become choppy, indicating those points when an 
accretion disk has formed (Figure |5| . 

These basic results were found to be relatively insen- 
sitive to the Mach number of the fiow or the presence of 
a small asymmetry in the upstream fiow. The original 
discovery of the fiip-fiop instability came from simula- 
tions in which a density gradient was imposed on the 
ambient medium. To see if the instability is affected by 
small gradients, we ran a simulation (model C) with a 
transverse density gradient as in equation ®, but with 
(r sin (j)) replaced by C as given in equation (pi . We used 
an asymmetry value of e = 0.06. We also ran a simu- 
lation (model B) with no imposed asymmetry but with 
a higher Mach number of 10. The results for these runs 
are shown in Figure [6] in the form of the evolution of 
the accreted specific angular momentum. In each case 
the growth of the fiip-fiop instability was similar to the 
canonical run with e = and a Mach number of 4. The 
growth rate of the instability in the higher Mach num- 
ber fiow was roughly 25% higher, but otherwise the three 
simulations were very similar. 

For applications to accreting neutron stars and black 
holes, we are interested in values of Rs that are far 
smaller than our canonical value of 0.037 i?^. Most sim- 
ulations in the literature are limited to relatively large 
accretors because decreasing Rg becomes computation- 
ally very expensive. The timestep in these simulations is 
set by the Courant condition near the surface of the ac- 




Fig. 4. — The flip-flop instability is shown in this time series 
of the gas density from our standard model. The early evolution 
seen in the first two frames is characterized by the swinging of the 
accretion shock from side to side, with ever increasing amplitude. 
Eventually the amplitude grows to the point where the accretion 
shock flips over the leading edge of the accretor, as seen in the last 
frame. During these phases when the accretion shock has wrapped 
around the accretor, the gas interior to the accretion shock forms 
a temporary accretion disk with nearly Keplerian flow. 



cretor, and hence scales as R, 



,3/2 



Some of our extended 
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runs with small accretors ran for more than 10^ time 
steps. In addition to our canonical model we present re- 
sults for models with an accretor radius of 0.0125 Ra and 
0.0037 i?a in Figure [7| 

A smaller accretor is more unstable. The growth rate 
is larger and the flow forms an accretion disk quicker. In 
the intermediate case the growth rate was roughly 50% 
higher and the first disk was formed half an oscillation 
earlier. In the smallest accretor case the growth rate was 
more than twice as high as the canonical run and the first 
disk formed in half the time - after only one oscillation 
of the accretion shock. The accretion disk also persists 
longer for smaller accretors. In the canonical model a 
temporary disk is formed for only brief epochs (typically 
one fiow time across the accretion diameter) with little 
impact on the quasi-periodic behavior of j{t). In the case 
of the smallest accretor a given disk phase can last for 
several fiow times, with very abrupt transitions between 
counter-rotating disks. 

The structure of the temporary accretion disks as a 
function of accretor size is illustrated in Figure [8] One 
can see a general trend of decoupling between the spiral 
shocks of the accretion disk and the leading bowshock 
as the size of the accretor is reduced. In model A with 
the largest accretor the leading bowshock is a direct ex- 
tension of a spiral shock reaching down to the surface of 
the accretor. In the intermediate model the spiral shock 
wraps at least one full revolution around the accretor be- 
fore it becomes a part of the external bowshock. In the 
case of the smallest accretor (model E) , the spiral shocks 
can be traced out to the bowshock, but it involves several 
windings. The result in this last case is a disk-like accre- 
tion fiow that is only weakly coupled to the bowshock. 

The oscillation period of the fiip-fiop instability in the 
small amplitude regime (before the formation of any ac- 
cretion disk) is relatively insensitive to the accretor ra- 
dius, with values comparable to the Keplerian period at 
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Fig. 5. — The growth of the flip-flop instabiUty from an initially 
equilibrium state is illustrated by the evolution of the specific an- 
gular momentum, j, of the accreted gas. The dashed curve in the 
top panel is the analytic fit described in equation (9). The mass 
accretion rate, shown in the bottom panel, is relatively constant 
except during brief periods of time when a disk forms and the mass 
accretion rate drops substantially. This low accretion phase is al- 
ways followed by a large jump, when the mass in such a disk is 
'flushed' onto the accretor resulting in a momentary spike in mass 
accretion. 
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Fig. 6. — The evolution of j is shown for models A, B (Mach = 
10), and C (e = 0.06). 
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Fig. 7. — The evolution of j is shown for an accretor of radius 
Rs = 0.0125 Ra in the top plot, along with j for the canonical 
model with Rs = 0.037 Ra (dashed line). The bottom plot shows j 
for an accretor of radius Rs = 0.0037 i?a, again compared against 
the canonical model. 

the accretion radius. The dependence of uji on Rs is hard 
to quantify because tOi decreases somewhat as the ampli- 
tude of oscillation increases. As a result, the value of uJi 
derived from the simulation data depends on the time 
window within which one fits the analytic function (eq. 
[9]). This window is forced to be early in the evolution 
for the small accretor case because of the rapid growth, 
but a better fit is found if a later window is used in the 
runs with a larger accretor. One can say, however, that 
the differences in uJi for different values of Rg are smaller 
than the variations in uji with oscillation amplitude. 

The mass accretion during the disk phase is driven by 
spiral shocks (Blondin 2000). This is illustrated in Fig- 
ure [9] where we show a streamline computed from the 
gas velocity at a fixed moment in time. Note that this 
streamline is not the same as the path taken by a fiuid 
element over time, but it provides a better illustration 
of the effects of spiral shocks on the flow. When or- 
biting gas encounters a spiral shock, its orbital velocity 
is decreased, resulting in the kinks seen in the stream- 
line. This process robs the gas of angular momentum. 
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TABLE 1 

Simulation Parameters 






Fig. 8. — Disk structure for three different sizes of the accretor, 
models A, D, and E. The vertical size of each image is Ra- The 
density scale is different in each of the three images to allow similar 
detail to be seen in each. 



Model 


Rs/Ra 


Mach 


7 


e 


UJr 


A 


0.037 


4 


1.333 





0.070 


B 


0.037 


10 


1.333 





0.088 


C 


0.037 


4 


1.333 


0.06 


0.070 


D 


0.0125 


4 


1.333 





0.104 


E 


0.0037 


4 


1.333 





0.170 


F 


0.037 


4 


1.667 





0.000 


G 


0.037 


4 


1.400 





0.068 


H 


0.037 


4 


1.500 





0.040 


I 


0.037 


4 


1.550 





0.031 


J 


0.0125 


4 


1.600 





0.011 


K 


0.00185 


4 


1.667 





0.008 



which is transferred to gas at larger radii through the 
pressure gradient behind the shock. The loss of angu- 
lar momentum in turn leads to an inward radial drift, 
and eventually accretion onto the central object. The 
specific angular momentum of the gas accreted during 
a disk phase is j ^ 1.18, independent of the size of the 
accretor. This is slightly larger than the Keplerian value 
of unity at the surface of the accretor, consistent with 
the eccentricity of the orbiting gas near the accretor. 




Fig. 9. — A streamline from a snapshot of the flow in the sim- 
ulation with Rs = 0.0037 -Ra illustrates the elliptical nature of the 
orbits and the influence of spiral shocks on the accretion of the or- 
biting gas. The shading represents gas pressure. The spiral shocks 
are marked by discontinuities in pressure as well as slight kinks in 
the computed streamline. 



Finally, we found that the adiabatic index has a dra- 
matic influence on this instability. We ran several simula- 
tions varying Mach number and density asymmetry with 
7=5/3 and Rs — 0.037 i?a, and in all cases the oscilla- 
tion of the accretion shock remained bounded. Typically 
the shock cone begins oscillation very early on but the 
oscillations never grow beyond a relatively small ampli- 
tude and no accretion disk is ever formed. This depen- 
dence on 7 is consistent with the results of Zarinclli et| 



al. (19951. They ran a simulation with a homogeneous 
upstream flow and 7 = 5/3 and reported steady flow for 
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an absorb ing i nner boundary. They also ran a simulation 
similar to FT with 7 = 4/3 and e — 0.005 and found the 
flip-flop instability, regardless of the inner boundary con- 
dition. The authors attributed this difference to e, but 
given the results reported here this difference is more 
likely due to 7. 
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1.55 



1.6 



1.65 



Fig. 10. — The growth rate of the flip-flop instabiUty is shown as 
a function of the adiabatic index of the gas, 7, for Rs = 0.037 Ra- 
The simulation with 7 = 1.60 exhibited small oscillations at a 
constant amplitude, but no evidence for growth. 



To further investigate this dependence on 7 we ran a 
series of simulations in which 7 was the only parame- 
ter changed. To quantify the instability we fit j(t) as 
in Figure [5] The results, shown in Figure |10[ show a 
monotonic decrease in the growth rate with increasing 
7. By 7 = 1.6 the accretion wake is no longer unstable. 
The wake still swings back and forth, but the amplitude 
of those oscillations is only a few degrees and it remains 
remarkably constant. 




with 7 = 1.6 for smaller accretors (see Table 1) and found 
the flow to be unstable, with increasing growth rate for 
decreasing Rg as in the case with 7 — 4/3. We also 
ran simulations with 7 — 5/3, searching to see if such 
flow is unstable if the size of the accretor is sufficiently 
small. With — 0.0025 i?^, the amphtude of the the 
oscillating accretion shock grew very slowly with cUr — 
0.006, but then seemed to level off before the amplitude 
became big enough for the shock to flip over the leading 
edge and produce a disk. The growth rate was about the 
same for a slightly smaller accretor, Rg — 0.00185 i?a, 
but as shown in Figure [iT] the amplitude was big enough 
to produce a flip after only six oscillations. 

4. CONCLUSIONS 

We have presented high-fidelity simulations of 2D pla- 
nar accretion onto a compact object in order to better 
understand the flip-fiop instability originally discovered 
by FT; From these simulations we can draw the following 
conclusions: 

1) First and foremost, we have shown that the flip- 
flop is a true overstability, where a small displacement of 
an equilibrium accretion shock leads to oscillations with 
exponentially growing amplitude. 

2) This instability does not require an asymmetry in 
the upstream flow; a result flrst suggested by [Matsuda] 
eFaLj^llQQTI . 

3) The growth rate is only weakly dependent on the 
Mach number of the upstream flow, with higher Mach 
number producing slightly faster growth of the instabil- 
ity. 

4) The growth rate increases with decreasing radius of 
the accretor. 

5) For a given accretor size, the stability of the flow 
is dependent on the adiabatic index of the gas, with the 
growth rate decreasing as 7 increases from 4/3. Planar 
accretion is stable for 7 > 1.6 in the case of a large 
accretor {Rs = 0.037 Ra), while 7 = 5/3 is only unstable 
for the smallest accretor used in our simulations, Rs — 
0.00185 



1 - 



0.5 - 




-1.5 I ' 1 ' ' 1 

20 40 60 80 100 

time (R^/V.) 

Fig. 11. — The growth of the flip-flop instability for an accretor 
of radius Rs = 0.00185 i?a and adiabatic index 7 = 5/3. 
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This critical value of 7, above which the flow is (at least 
marginally) stable, depends on the size of the accretor. 
As shown above, a smaller accretor produces a faster 
growth rate of the instability. A smaller accretor also 
means a larger critical value of 7. We ran simulations 
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